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ABSTRACT 

We present results concerning adiabatic inertial-mode oscillations of non-stratified 
superfiuid neutron stars in Newtonian gravity, using the anelastic and slow-rotation 
approximations. We consider a simple two-fluid model of a superfiuid neutron star, 
where one fluid consists of the superfiuid neutrons and the second fluid contains all the 
comoving constituents (protons, electrons). The two fluids are assumed to be "free" 
in the sense that vortex- mediated forces like mutual friction or pinning are absent, 
but they can be coupled by the equation of state, in particular by entrainment. The 
stationary background consists of the two fluids rotating uniformly around the same 
axis with potentially different rotation rates. We study the special cases of co-rotating 
backgrounds, vanishing entrainment, and the purely toroidal r-modes, analytically. 
We calculate numerically the eigenfunctions and frequencies of inertial modes in the 
general case of non co-rotating backgrounds, and study their dependence on the rel- 
ative rotation rate and entrainment. In these non-stratified models we find avoided 
crossings only between associated mode-pairs, eg. an "ordinary" mode and its "su- 
perfiuid" counterpart, while other mode- frequencies generally cross as the background 
parameters are varied. We confirm (for the first time in a mode calculation) the onset 
of a "two-stream instability" at a critical relative background rotation rate, and study 
some of the properties of this instability for the inertial modes. 
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>^ ■ 1 INTRODUCTION 



The oscillations of rotating compact stars is a subject that 
has attracted interest for a considerable time. This is natural 
since the associated issues range from fundamental applied 
mathematics (eg. the stability of rotating self-gravitating 
fluid configurations) , to mainstream astrophysics (eg. Helio- 
seismology and attempts to infer the Sun's rotation profile 
from observed modes of oscillation) and exotic neutron star 
physics (eg. the gravitational-wave driven instability of the 
r-modes and various viscous damping mechanisms, like hy- 
peron bulk viscosity). To date, most investigations have as- 
sumed that a rotating star can be appropriately described 
by a perfect fluid model. While such models are relevant in 
many contexts, they do not provide an adequate description 
of mature neutron stars. Once a neutron star has cooled be- 
low 10^ — 10^" K, i.e. within minutes to months after its 
birth, it's outer layers will form a crystalline lattice of nu- 
clei. At the same time, the fluid core is expected to contain 
several superfluid/superconducting components. 

This paper concerns the dynamics of rotating superfiuid 
neutron stars. In particular, we study the inertial modes of 
a simple two-fluid model appropriate for the conditions that 



prevail in the outer core of a neutron star. The two fluids, 
which represent superfiuid neutrons and a congomerate of 
all comoving constituents (protons, electrons), are coupled 
via the equation of state (in particular via entrainment), 
but are otherwise allowed to move at independent velocities. 
Our background model is a stationary two-fluid configura- 
tion with constant entrainment, with the two fluids rotating 
uniformly around the same axis with rotation rates and 

Previous studies of superfiuid inertial modes (includ- 
ing preliminary studies of the zero-frequency subspace), 
namely Lindblom & Mendell (2000); Sedrakian & Wasser- 
man (2000); Andersson & Comer (2001); Comer (2002); Lee 
& Yoshida (2003); Yoshida & Lee (2003a,b), have aU been 
restricted to co-rotating backgrounds Sin = Sip. This study 
is the first to allow for the general case of a background with 
two fluids rotating at different rates. This is expected to be 
the quasi-stationary "ground state" of a superfluid neutron 
star due to its emission-induced spindown and the weak cou- 
pling to the superfluid components. This background model 
is also the starting point of all viable models of Vela-sized 
glitches. In other words, by allowing for difi'erent rates of 
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rotation, wc have taken a crucial step towaxds more realistic 
modelling of the dynamics of mature neutron stars. 



2 THE TWO-FLUID NEUTRON STAR MODEL 

We take as our starting point the "standard" two-fluid 
model for superfluid neutron stars (e.g. Lindblom & Mendell 

(1994); Lcc (1995); Langlois ot al. (1998); Prix & Rioutord 

(2002) ), in which one assumes that the protons and elec- 
trons are locked together by the magnetic field and viscosity, 
while the superfluid neutrons form an independent fluid due 
to their lack of viscosity. Our model neglects the presence of 
the clastic crust as well as the potential presence of exotic 
matter in the deep neutron star core. In essence, the model 
is expected to be relevant for the outer neutron star core. 
By studying the global modes of oscillation of this model, 
we hope to gain insight into the complex dynamics of any 
two-fluid system. Even though we are not considering a de- 
tailed realistic neutron star model (the construction of which 
would be very difficult given our current level of understand- 
ing) we expect to learn much about qualitative aspects that 
should remain relevant also in more complicated settings. It 
is also interesting to note, cf. comments made by Sedrakian 
& Wasserman (2000), that the study of two- fluid models 
may be of significance in laboratory contexts, for example 
in the study of rotating heavy nuclei using the compress- 
ible liquid approximation of the Bohr- Wheeler model, or for 
rotating mixtures of Bose- Einstein condensates. 

A general Newtonian formalism to describe mixtures of 
charged and uncharged fluids has been developed by Prix 
(2003a, b), based on a variational principle that was first 
developed in a fully relativistic framework by Carter and 
coworkers (Carter 1989; Carter & Khalatnikov 1992; Carter 
& Langlois 1998). In particular, Prix (200,3a) developed a 
general two-fluid neutron star model allowing for tempera- 
ture gradients and dissipation through mutual friction and 
/3-reactions between the two fluids. For the present appli- 
cation, however, we assume a "cold" neutron star in which 
we can neglect temperature effects, so we set T = Q, and 
we also neglect mutual friction and non-adiabatic processes 
like /3-reactions. The resulting framework, which is identi- 
cal to that used by Prix & Rieutord (2002) and Andersson 
& Comer (2001), is briefly introduced in this section. Note 
that although the formalism used here is different from the 
one more commonly found in the Newtonian literature (e.g. 
Lindblom & MendeU (1994); Lee (1995); Lee & Yoshida 

(2003) ), which is based on the "orthodox" superfluid for- 
malism introduced by Landau, the two frameworks can be 
shown to be strictly equivalent, as discussed by Prix (2003a). 

Our two-fluid model consists of a neutron and a "pro- 
ton" fluid (the latter actually consists of the comoving pro- 
tons and the electrons). Therefore the kinematic variables 
are the particle number densities rin and rip together with 
the respective transport velocities Wn and Up. The corre- 
sponding transport currents are naturally expressed as 



nx = nx vx , 



(1) 



where X G {n, p} is the constituent index (the repetition of 
which does not imply summation). An important quantity 
for our analysis is the relative velocity A between the two 
fluids, which we deflne as 



A S Vp - Vn . (2) 

The dynamics is governed by the Lagrangian density 

(3) 



1 2 1 2 



2 " 2 

where p = mnTin + rUpUp is the total mass density, 3> is the 

gravitational potential and £ the energy function or "equa- 
tion of state" of the system. The general form of the equation 
of state is £ = £{nn, rip, A^), which determines the first law 
of thermodynamics in the form 

d£ = /i" dun + drip + a dA.^ , (4) 

defining the chemical potentials /n" and iJ', as well as the 
entrainment a. The conjugate momenta for the two fluids 
are defined by the total differential of the Lagrangian density 
A, namely 

dA= ^ [p^ -dnx + ipo -m^^)dnx]-pd^, (5) 

X=n,p 

In the following we assume the two masses to be equal, so 
we set mP = m" = mb. With the explicit form (3) of the 
Lagrangian and the first law (4), we can express these con- 
jugate momenta as 



m-b {Vn -I- £n A) , 



(6) 
(7) 

(8) 
(9) 



p = 

p'" = m-b {vp - £p A) , 

n n I 1 2 n 

Po = -A* + ^m^v^-Vn-p , 
Po = -p.^ + ^rnbvl-Vp-p^ , 

where we have defined the dimensionless parameters ex 
characterizing entrainment by 

£x = ■ (10) 

rUbUx 

Sometimes it is more convenient to use a single entrainment- 
parameter e, which we choose to be Sp, so we have 



Ep = e , and Sn ■ 



1 



(11) 



in terms of the proton fraction Xp, which is naturally defined 
as 

Xp = — , with n = rin + np . (12) 
n 

We note that this definition of the entrainment e is different 
from another definition, e, which is sometimes found in the 
literature (e.g. Lee & Yoshida (2003); Lindblom & Mendell 
(2000)). The relation between these two different definitions 
is simply (see Prix et al. (2002) for further discussion) 

6 = ^''^ . (13) 
Un — en 

We assume that the timescale of oscillations is much shorter 
than that of /3-reactions. Therefore strict conservation of 
neutrons and protons applies, i.e. we have 

0, 
0. 



dtUn + V • (rinVn) 

dtUp + V • {ripVp) 



(14) 
(15) 



As shown by Prix (2003a), the equations of motion for the 
two fiuids can be derived from the Lagrangian density (3) 
using a "convective" variational principle. They can be writ- 
ten in the form 
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{dt +VX-V 



P +Pi Vi-x 



11 

nx 



(16) 



where is the "external" force density acting on the fluid 
X, and the scalars Q'^ are defined as 

QX X X jxL I X 

= po - m 9 + vx ■ p 



= -1^^ + ^rrinVx - rrin^ . 



(17) 



In the absence of "external" forces acting on the whole sys- 
tem, the hydrodynamic force densities in (16) have to 
satisfy /" + Z'' = as a Noether identity of the variational 
principle. This still allows one to describe a mutual force 
/mut acting between the two fluids. It could be caused, for 
example, by collisions of the electrons with the neutron vor- 
tices (e.g. cf. Alpar et al. (1984)). Such a model would be 
characterized by /° = — = /mut- As a first step, however, 
we only consider the "free" limit and postpone the inclusion 
of mutual friction and viscosity to future work. Our "free" 
model is therefore characterized by 



f"" =0, 



for X = n, p . 



(18) 



3 THE STATIONARY BACKGROUND 

We assume the background to be stationary and axisymmet- 
ric, with both fluids rotating around the a-axis with rotation 
rates Qn and fip respectively. Hence 

vx = ^ix f , and A = (f2p — Qn) V , (19) 
where (p is the axial Killing vector, given by 



(20) 



In spherical coordinates, i.e. a;* £ {r, d, (p}, this vector 
has the components y?' = (0, 0, 1), and its norm is 

(p^ipi — s\v? 6. With the cntrainmcnt relations (6) and (7) 
we can now write the background fluid momenta as 



p" = mb (fin - £n(^^n - fip)) , 

pP = rrih (fip — ep(np — fin)) V ■ 



(21) 
(22) 



In the following it will be convenient to introduce as a short- 
hand notation the tilde-operator acting on a constituent 
quantity, Q.x say, as follows 

hx=^x-ex{0.x-nY), where F 5^ X . (23) 

This allows us to rewrite the background-momenta as 

= ruh^x <p ■ (24) 

We restrict our attention to models with uniform rotation, 
i.e. Vilx = 0, therefore the background vorticitics are 

x 



V X p = 2miyQ,x z + (fix — fiy) x Vex 



(25) 



where z is the unit vector along the z-axis. We see that in 
the general case of a varying entrainment ex and different 
background rotation rates f2n fip, the vorticities are no 
longer aligned with the rotation axis. In other words, they 
acquire a non-zero ^-component and the system is in a state 
which resembles differential rotation. 

As a first step towards a complete understanding of the 
dynamics of rotating two-fluid systems, we will focus on one 
of the simplest possibilities. We make the assumption that 



the entrainment ex is constant throughout the star, which 
means that we have Vex = and therefore also Vfix = 0. 
We note that assuming both entrainment parameters ex to 
be constant also requires a constant proton fraction Xp. In 
the following we therefore consider a non- stratified neutron 
star model (i.e. Vxp = 0) with constant entrainment. This 
model is admittedly simplistic, but as we will see in the 
following, it nevertheless allows for a rich phenomenology. 



4 LINEAR OSCILLATIONS 

4.1 Oscillation equations in harmonic basis 

Assuming uniform rotation and a constant entrainment 

model as discussed above, the linear perturbation of the 
equations of motion (16) can be obtained in the form 

{dt + Q.x£^)— + 2Q.xCx+Vil)x=Q, (26) 
mb 

where the Lie derivative explicitly gives 

£if, 5pi = ip^'^ j5pi + 5pj'^ np\ and 5 represents an Eu- 
lerian perturbation. We have also defined the Coriolis-term 
Cx for each of the two fluids as 

Cx = zx5vx, (27) 

and a scalar potential i/'x representing the "effective" pres- 
sure perturbation, namely 

ruhipx = Sn^ +mh5(t>+ (p'^ — mbVx) • 5vx ■ (28) 

The background is assumed to be stationary and axisym- 

mctric, so we can look for eigenmode solutions of the form 
^i(;^t+mcp)^ and the comoving time derivative in (26) can be 
directly replaced by 



{dt + (lx£^) i{ui + mClx) ■ 



(29) 



The practical advantage of using the Lie-derivative £^ here 
is that the substitution £^ im holds for any geometric 
object (e.g. a vector as in Eq. (26)) with a i/j-dependence of 
the form e'""'', while this is only true for the simple direc- 
tional derivative if ■ V when it is applied to scalars. 

Linear perturbation of the conservation equations (14) 
and (15) leads to 



dtSnx + V ■ {nxSvx + Snxvx) = . 



(30) 



In the present analysis we are only interested in inertial 
modes, which are characterized by frequencies of the or- 
der of the rotation rate fi. Since this frequency is usually 
much lower than that of the lowest order p-mode frequency 
iUp, we can simplify the problem by using the anclastic ap- 
proximation (which effectively "filters out" the p- modes). 
As discussed in more detail in Appendix A, the anelastic 
approximation consists of replacing the conservation equa- 
tions by 



(31) 



V • (nx Svx) = + [ — 



The lowest-order p-modc frequency cOp is of the order of the 
sound-crossing frequency, i.e. ujp = 0{cq/R), where Co is an 
averaged sound speed and R is the stellar radius. For low- 
frequency modes such as the inertial modes we can therefore 
drop the higher order corrections in (31), which account for 



© 0000 RAS, MNRAS 000, 000-000 



4 R. Prix, G.L. Comer and N. Andersson 



the "elasticity" (i.e. compressibility) of matter. In the follow- 
ing we also restrict ourselves to slowly rotating backgrounds. 
Since 



(32) 



the star remains spherical if we neglect the centrifugal de- 
formation. It is important to note the difference between 
the slow-rotation approximation, which compares the rota- 
tion rate Q, to the Kepler limit Q.k = 0(%/47rGpo), and the 
anelastic approximation, which is relevant for mode frequen- 
cies that are small compared to Up. 

In order to solve the perturbation equations (26) and 
(31) for inertial modes, we first transform them into a pseudo 
one dimensional problem by expressing all angular depen- 
dencies in terms of the the spherical harmonics Yi^{9,ip). 
The spherical harmonics are the eigenvectors of the angular 
Laplacian, namely 



(33) 



Since these functions form a complete orthonormal basis, we 
can expand a scalar field tpx as 

iJx{r,e,y,)=ij'x{r)Yr{9,v), (34) 

where here and in the following automatic summation over 
repeated "angular" indices {I, m,...) applies. Because of the 
assumption of an axisymmetric background, the various m- 
contributions can be decoupled, and therefore we can con- 
sider each value of m separately. In order to express a vec- 
tor field in a similar manner, we use the "harmonic basis" 
{Ji™, S™, T["}, which is defined in terms of the spherical 
harmonics as 

RT = Vr , ST = Yr , TT = ST X Vr . (35) 
We can then expand the velocity perturbations 6vx as 



5v 



■RT + vjcir) sT-mrr 



(36) 



Using the entrainment relations (6) and (7), and assuming 
a constant entrainment model we obtain 

= dvx - ex {dvx - Svy) = Svx , (37) 

nib 

where we have used the definition (23) of the tilde-operator. 
This can conveniently be written in the harmonic basis as 



TOb 



■RT + v^ sr 



(38) 



The gradient of a scalar field (34) is readily expressed as 

(39) 



V^x^^k {r)RT+^Ur)Sl' 



where the prime represents a radial derivative. The expres- 
sion for the Coriolis-terms (27) in the harmonic basis is 
found after a straightforward, but somewhat laborious, cal- 
culation to be 



-I (mV^ + Plu'x) 



(40) 



where we sum over the repeated "angular" indices (fc and 
/), and the constant matrices /3f and 7^ are defined as 

/3f = lQi+i6k,i+i- il + l)QiSk,i-i, (41) 

-yi = (e -l)QiSk,i-i + lil + 2)Qi+r6k,i+i, (42) 

with the usual definition 



Q 



4/2 _ 1 



(43) 



4.2 The general eigenmode equations 

Putting all the pieces together, we can now express the com- 
plete system of equations (26) and (31) in the harmonic basis 



KxVx 



KxUy, 



1 + r^) W^i 

nx 



i{i + im 



1(1 + 1) 



KxWk - 2{mVk + Aui) = 
(mW'x+mVjc+jiU^) = 



0, 

2rVry 



/(/ + 1) 
where we defined 
u> + mQx 

Kx ~ 



ipl^W^+mU'x+'ylV^) = 0, 



and 



■'tpx ■ 



(44) 

(45) 
(46) 

(47) 



(48) 



nx 2nx 

Wc note that this definition of the dimcnsionlcss fre- 
quencies nx reduces to the usual single-fluid deflnition 
K= [uj + mfi) /Q. in the case of comoving fluids, or in the 
absence of entrainment. In both of these cases we have 
Q,x — > fix as seen from the definition (23). 

The boundary conditions at the centre of the star 
(r = 0) consist of the regularity requirement of the harmonic 
expansion (34) and (36), which implies the asymptotic con- 
ditions 



W'x 



W'x 



0{r) as r^O. 



(49) 



At the surface (r = R) we require another regularity condi- 
tion due to the divergent term n'x/nx in the conservation 
equations (44). As discussed in the Appendix A, this is a 
consequence of the anelastic approximation. The resulting 
surface boundary condition is therefore 



Wx{R)=0, 



(50) 



i.e. the radial displacement vanishes at the surface. For mod- 
els with a vanishing surface density we do not need to impose 
an explicit condition on the pressure variation at the surface, 

as we have 



AP = 6P + ^ ■ T7P = SP = ^nx Sn^ 



(51) 



The vanishing of the Lagrangian pressure perturbation AP 
is therefore ensured provided that the Sfi'^ arc regular at the 
surface. As our rmmerical scheme can only find such regular 
solutions, this boundary condition is implicitly guaranteed 
to hold. 



4.3 Special case: zero entrainment 

We note that the only coupling between the neutrons 
[X = n) and the protons {X — p) in the eigenvalue system 
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(44)-(47) is caused by the entrainmcnt, cf. the definition of 
the tilde-operator (23). In the case of zero entrainment, i.e. 
£n = Ep = 0, we obtain two uncoupled eigenvalue systems. 
Both systems arc formally identical, and therefore both have 
the same solutions for kx, i.e. 



(52) 



where Kord represents the single-fluid solutions. However, 
from the definition (48) we see that these correspond to dif- 
ferent mode frequencies when the rotation rates of the two 
fluids are different, i.e. 



Wn = (/tord — m) On , and 



(Kord -m) dp, (53) 



which implies that these two solutions cannot form a single 

mode solution when Qn 7^ f2p. The two modes in this case 
therefore correspond to only one of the two fluids oscillating 
while the other fluid is at rest, i.e. 



u> = u)p : Sv„ = , 



and 5vp = , 
and 5vp . 



(54) 
(55) 



From the fact that one of the two fluid amplitudes neces- 
sarily vanishes when e — > we deduce that the correspond- 
ing amplitude will actually change sign at this point. Wc 
can therefore conjecture that if the two fluids were predom- 
inantly in phase before crossing e = 0, then they will be 
predominantly in counter-phase afterwards and vice-versa. 
We will see in Sect. 6 that our numerical results agree per- 
fectly with this conjecture. 



4.4 The r-mode sub-class 

The subclass of purely cixial inertial modes, commonly re- 
ferred to as r-modes, has generated a lot of interest due 
to its strong instability with respect to gravitational waves 
(Andcrsson 1998; Friedman & Morsink 1998). Therefore it 
is interesting to see if this subclass still exists in the su- 
perfluid case, and how its properties are modified. A purely 
axial velocity perturbation is proportional to TJ" , so we set 
W'x = Vk = 0. In this case the equations (44)-(47) reduce 
to 



2nxAU^+irij'x = 0, 

2nx'riu^+iiii + i)i^'x = 0, 

0. 



^Kxl{l + l)U'x-mUlc 



(56) 
(57) 

(58) 



While the first two equations (56) and (57) allow one to cal- 
culate the eigenfunctions Ux and rpx > they do not constrain 
the eigenvalue in any way. The third equation (58), however, 
leads to an algebraic constraint for the existence of a non- 
trivial axial solution. In order to find this constraint, we use 
the explicit expressions for (58) for the two fluids, i.e. 



(59) 
(60) 



^l{l + l)Kn[{l-er,)U^ + enUl]-mUL = 0, 

^lil + l)Kp[{l-Sp)Ul + epUL]-mUl = 0, 

from which we can eliminate the Ux (assumed nonzero) 
to obtain the following dispersion relation for superfluid r- 
modes 

1(1 -h 1)(1 - en)(w + mf2n) - 2mfin 



X l{l -h 1)(1 - £p){u) + mQp) - 2mnp 

— 1^(1 + l)^£nep{u! + miln){LJ + mQp) = . 



(61) 



We note that this corrects the dispersion relation that was 

used by Andcrsson ct al. (2003) in a discussion of the super- 
fluid two-stream instability of the r-modes. 



5 THE CO-ROTATING CASE On = Qp 

Before turning to the rmmerical solution of the general case 
with fin 7^ f2p, it is instructive to study the special case 
of two co-rotating fluids, where wc have fix = O. The lin- 
earized perturbation equations (26) and (31) then take the 
form 



in [Svx - £x{Svx - Sw)] + 2Cx + 



V^x 



= 0, 



V-{nxSvx) = 0, 



(62) 
(63) 

(64) 



where we have defined 
u! + mil 

which is the usual dimensionless frequency of inertial modes 
in the co-rotating frame. It is interesting to see under which 
conditions this system can be separated into purely co- and 

counter-moving modes. We therefore introduce the usual 
variables corresponding to these two mode-classes, namely 



5A: 

Sv ■ 



SVp - 6Vn , 

Xp SVp + Xn SV„ 



6p- 
5/i ; 



/ (65) 

Xpipp + XntPn ■ 



In terms of these variables the oscillation equations can be 

rewritten as 



V- {nSv) 
IK 5v + 2z X 5v + Vip 
V • {nXnXp SA) 
i-y'^ k5A + 2z X SA + V6/3 

where we have defined 

1 

7 = 



0, 

6(3 Vxp , 
—nSv ■ Vxp , 
0, 



1 £n £p 



(66) 
(67) 
(68) 
(69) 



(70) 



We see that the variables {5t;, 5/^}, which are characteristic 
of "ordinary" -type modes, decouple from the "superfluid" 
variables {5A,(5/9}, if and only if the background model is 
not stratified, i.e. if Vxp = 0. This is exactly the same con- 
dition that was found in the case of a static superfluid neu- 
tron star (Prix & Rieutord 2002). We further see that in the 
non-stratified case the equations governing the two mode- 
families are equivalent, and the "ordinary" -type mode fre- 
quencies /tord are therefore related to the "superfluid" -type 
ones by 



^ 7 '^ord , 



(71) 



It is well known (Bryan 1889) that the inertial mode frequen- 
cies K f2 of an incompressible fluid are bounded (and form a 
dense set) in the interval [—20, 20] . In the compressible case 
we still expect this to hold approximately. This will therefore 
also be true for the "ordinary" -type modes in the co-rotating 
case, but relation (71) shows that the corresponding inter- 
val for the "superfluid" -type modes will be governed by the 
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factor 7. This scale factor depends only on the proton frac- 
tion Xp and the entrainment e, and can in principle take 
any value between [—00, +00]. For e < we have 7 e (0, 1), 
i.e. the "suporfluid"-type mode frequencies lie closer to the 
origin than their "ordinary" -type counterpart, and they are 
bounded by a smaller interval than the "ordinary" modes. 
For e > on the other hand, the "superfluid"-type mode 
frequencies lie further away from the origin than their "ordi- 
nary" counterparts and their bounding interval is larger. Wc 
also note that the "ordinary" -type modes are independent of 
the entrainment, as expected from their strictly co- moving 
character. If we express the scale factor (70) in terms of the 
alternative entrainment parameter e as defined in (13), we 
find 



7 = 1 + 



(72) 



Therefore we see that the "superfluid"-type mode frequen- 
cies are linear in e. This has been found previously for the 
r-mode subclass by Andersson & Comer (2001). It was also 
observed numerically for inertial modes by Lcc & Yoshida 
(2003) and Yoshida & Lee (2003a), although only as an "al- 
most" linear dependence. This slight discrepancy is not sur- 
prising as their background model is stratified (i.e. Xp is not 
constant) and the above decoupling of the mode-families is 
therefore expected to hold only approximately. 



Substituting this eigenfunction and the eigenvalues (75) 
back into (58), we find the following relation between the 
two amplitudes 



K — Kord 
K = 7 ^ord 



Cp — Cn , 

n — 



(80) 
(81) 



which corresponds to purely co-moving and counter-moving 
r-modes respectively. 



6 NUMERIAL RESULTS FOR THE GENERAL 

CASE On / Op 

In the following wc choose the proton rotation rate Op as 
the "reference" rotation rate. This choice is motivated by the 
fact that the observed rotation of neutron stars (via pulsar 
emission) is thought to be related to the charged components 
(assumed to co- rotate with the crust), while the rotation rate 
On of the superfluid neutrons is not directly observable. We 
define the relative rotation rate 7?. as 

^ _ On Op 



Op 

With these definitions we can write (23) as 

On = Op[l + (l-en)7^] , 

Op = Op [1 Ep 7e] . 



(82) 



(83) 
(84) 



5.1 The r-modes of the co-rotating model 

In the case of a co-rotating background, the r-mode disper- 
sion relation (61) reduces to 



[(1 - En) K - Kord] [(1 - £p) K - Kord] 

where 

2m 



/^ord 



/(/ + i) 



0, (73) 



(74) 



is the standard single-fluid r-mode frequency. Solving the 
quadratic dispersion relation we find that the two r-mode 
frequencies of the superfiuid problem are 



^ — {^ord) T'^ord} i 



(75) 



in agreement with the general relation (71). 

In order to find the corresponding eigenfunctions, we 
eliminate V'x from (56), (57) to obtain 



+ 



« + 2, 



l^ 



0. (76) 



This equation has to hold for every I > m, and we can there- 
fore extract the two simultaneous conditions 



Q 



(-1-1 



+ 1) ■ 



r 



= 0, 
= 0, 



(77) 
(78) 



for which the only non-trivial and non-singular solution is 



Cxr 



with m = I . 



(79) 



Further introducing 

1 



— — , and 2i/„ = — ■ , 

l + (l-en)7^ ^ l + SpTZ 



we can express nx defined in (48) as 



= 2Vp Ko , 



(85) 



(86) 
(87) 



where kq is a dimensionless "reference frequency" of the 
mode, which we define as 



a; -|- mOp 



(88) 



in analogy with the usual single-fiuid definition. With these 
definitions we can write the system of equations (44)-(47) 
as a one dimensional infinite eigenvalue problem for kq in 
the form 



(89) 



where the Ai and Bi are linear operators and is the eigen- 
vector 

¥ = l^wL wl, k', vl, ul ul, v^i, Vi'p} . (90) 

The explicit form of these equations is given in Appendix B. 
By taking the sum over I only up to a finite value Zmax, 
we can solve the resulting finite eigenvalue problem using 
the LSB spectral solver, which is based on the efficient in- 
complete Arnoldi-Chebychcv algorithm. This is the same 
method that Prix & Rieutord (2002) used to study non- 
radial oscillations of non-rotating stars. 

Most of the numerical calculations in the following have 
been performed for both a uniform-density background (i.e. 
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polytropic index A'^ = 0) and a polytropic background with 
N = 1 (for each of the fluids). The results are quite similar 
and we therefore only present the polytropic case here. Fur- 
thermore, we only considered the case m = 2, which is ex- 
pected to be the most relevant for gravitational-wave emis- 
sion. The results for higher values of m are not expected 
to show any qualitative differences. In all of the following 
sections except for Sect. 6.4, we use a neutron star model 
with "canonical" values a;p = 0.1 for the proton fraction and 
e = 0.6 for the entrainment, which conveniently results in a 
scaling-factor (70) of 7 = 3. In Sect. 6.4 on the two-stream 
instability we choose these values to be Xp = 0.2 and e = —2, 
which leads to 7 « 0.2857. 

For the low-order inertial modes considered in this pa- 
per a radial resolution of 30 Chebychev polynomials and an 
angular resolution of about 10 spherical harmonics is used 
in most cases, which proves sufficient to obtain a numerical 
precision of the order of lO"'^. In the case of r-modes, we 
compared our numerical results to a direct evaluation of the 
dispersion relation (61) and found an agreement better than 
10"^ in all cases considered. In the co-rotating case {TZ = 0) 
our numerical results for the "ordinary" modes agree per- 
fectly (up to the given precision ~ 10~^) with the single-fluid 
results in the literature (e.g. Lockitch & Friedman 1999), and 
the "superfluid" modes satisfy the relation (71) as expected. 



6.1 Angular convergence and inertial-mode 
labelling 

It was shown by Lockitch & Friedman (1999) that for m 7^ 
the lowest non-zero Z-coefBcient in the harmonic expansions 
(34) and (36) is necessarily I = \m\. Furthermore, in the 
case of a uniform density background model (i.e. A'' = 0), 
it is known that the harmonic expansion (36) of the solu- 
tion 5v stops at a finite Iq. In fact, the corresponding coef- 
ficients can be calculated analytically. There are always ex- 
actly j = lo — \m\ -\- 1 mode-solutions for any given m 7^ 
and Zo > In the case of a polytropic background with 
N = 1 the solution turns out to be quite similar to the uni- 
form case, except that the expansion does not stop after a 
finite number of terms. Instead it converges exponentially 
beyond I = lo- This behaviour is illustrated in Fig. 1, which 
shows the angular expansion coefficients of the (axial-led, 
m = 2) inertial mode with k = —1.308 for N = and the 
analogous mode k = —1.43392 for A'^ = 1. We see that in the 
uniform background case there is a sharp drop after j = 5, 
as the higher-order coefRcients are analytically zero, while 
in the polytropic case we observe an exponential falloff. The 
quantitative and qualitative similarity to the uniform model 
allows one to associate the modes of the polytropic model 
with corresponding modes of the uniform model. In case of 
doubt it should always be possible to associate modes via 
a continuous transformation of A'^. We can therefore conve- 
niently label the modes by their "quantum numbers" m, j 
and an additional index n G [1, j] accounting for the j differ- 
ent solutions at given m and j. As a convention we choose 
to order the modes by increasing frequency n, so we label 
the inertial modes as 



< 



(91) 



o 



n I r 



N=0 
N=l 



I L_ 

7 8 9 10 

l-m + 1 



Figure 1. Angular convergence of the ("ordinary") inertial mode 
with m = 2 and j = lQ — m + l = 'o (which we label at 
7^ = for the homogeneous background model {N = 0) and the 
A'^ = 1 polytrope. The plotted quantity c; represents the magni- 
tude of the harmonic expansion coefficients (34) and (36) of the 
eigenmode. 
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where the inequalities obviously refer to the eigenfrequency 
of the corresponding mode. 



Figure 2. The polar-led "ordinary" inertial mode 2^°^^ for Xp = 
0.1 and e = 0.6 for a co-rotating background 7?. = (left panel) 
and 7?. = 0.1 (right panel). The next higher /-contributions are 
one order of magnitude smaller and are not included in this graph. 



6.2 The effect of relative rotation TZ 

We have seen in Sect. 5 that, in the co-rotating case the 
inertial modes of non-stratified stars can be separated into 
purely co- and counter-moving families. This is no longer 
true when we allow for a non-zero relative rotation 7?. 7^ 0. 
Similar to stratification (cf. Prix & Rieutord (2002)), the 
relative rotation introduces a coupling between these mode- 
families, leading to a deviation from the strictly co- and 
counter-moving nature of the modes. This is shown in Fig. 2. 

Compared to the effect of stratification, however, the 
mode-coupling induced by the relative rotation TZ (in the 
absence of stratification) seems to be of a much weaker 
nature. Although the two fluids are no longer strictly co- 
or counter-moving, they always have a well-defined phase- 
relation, in the sense that they are either strictly in phase or 
in counter-phase. Changing TZ does not change the position 
of the nodes of the mode. This can be seen in Fig. 3 which 
illustrates the transition of the 2l°2^ mode being in phase to 
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Figure 3. Neutron- and proton-amplitudes for different relative 
rotations TZ. The plot shows the I = 2 components Un and Up of 
the mode ^1°^^ for n = Q, -0.1, -0.4, -0.6. The normalization 
is such that Up{R) = 1, for which Up is seen to be invariant under 
changes of H. 
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n 

Figure 4. Frequencies kq of "ordinary" (— ) and "superfiuid" (- -) 
m = 2 inertial modes as functions of the relative rotation rate TZ 
for an N = 1 polytropic background model with Xp = 0.1 and 
£ = 0.6. The two common crossing points (94) are marked by 'x'. 
The modes presented here are the 6 lowest order inertial modes, 

2^(1) to 2^(3)- 



being in counter phase when varying TZ. Furthermore, this 
couphng does not lead to general avoided crossings between 
inertial modes, as can be seen in Fig. 4 in which we show 
the mode frequencies of the lowest-order inertial modes as 
functions of the relative rotation rate 71. 

A striking feature of this graph is that there are two 
common crossing-points for the mode frequencies. This can 
be understood as follows: the system of equations (B1)-(B8) 
has a singularity when z/x — > oo for X = n or X = p. We 
see from (85) that this happens at the relative rotation rates 
TZ — — j-i— and TZ = —j-, respectively. At these singular 
points the system of equations reduces to the following con- 
straints: 

Vn^oo : Ko = -mn , or Wl = - ^° Wl , (92) 

1 £n 



Vp^oo : Ko=0, or wl, = ---^^Wl, (93) 

1 Sp 

and similar amplitude constraints hold for Vx and Ux. The 
two common crossing points therefore have to be 

(7^, Ko) = {(--, 0), (-^,-m7^)}. (94) 

£p 1 £n 

These points are are marked by 'x' in our various frequency 
plots. The solutions at these critical relative rotation rates 
fall into two classes: modes that cross at the common cross- 
ing point, and modes that satisfy the amplitude relations 
(92) or (93). These analytical results agree perfectly well 
with the numerical findings and provide a good consistency 
check of our numerical results. We note that while in Fig. 4 
it seems as if each of the modes necessarily passes through 
one of the two crossing points, this is not generally the case, 
as will be seen in Fig. 7 for a different choice of parameters. 

It is interesting to note that the two critical relative 
rotation rates correspond to the vanishing of the angular 
momentum of one of the two fluids, i.e. f-n ^ oo corresponds 
to p" = and ^ oo is equivalent to = 0. This is 
obviously an eflfect of the entrainment: the fluid is rotating 
but has zero angular momentum! As a result, the Coriolis 
force acting on this fluid vanishes and the mode becomes 
stationary in the reference frame of the respective fluid. As 
we have chosen Q,p as our reference rotation, we find kq = 
for p'' = 0. The nonzero crossing frequency ko = —mTZ for 
p" = simply corresponds to the zero frequency in the 
neutron-frame observed in the proton-frame. 

While there are no general avoided crossings, the cou- 
pling induced by TZ does lead to avoided crossings between 
corresponding "mode-pairs", i.e. between the "ordinary" 
mode and its "superfluid" counterpart, as can be seen in 
Fig. 5. We note that the labelling Kord and Kgf used in 4 to 
refer to "ordinary" or "superfluid" modes is defined by con- 
tinuing the mode from TZ = 0. This labelling is somewhat 
arbitrary, however, as for 7?. it does not reflect the co- 
or counter-moving nature of the mode. Neither does it imply 
the mode to be in phase or in counter-phase, as can be seen 
from Fig. 3. In Fig, 5 and in the following it will often be 
more interesting to indicate the phase-character of a mode, 
so we will write k+ for modes with in-phase fluid motion, 
and K_ for modes where the fluids are in counter-phase. 

As we have already seen in Fig. 3, the relative phase is 
not an invariant property of the "ordinary" or "superfluid" 
mode families. For example, in Fig. 5 the "superfluid" modes 
are always in counter-phase, while the "ordinary" mode is 
in phase in a certain region but in counter-phase in another. 
We note, however, that the ordinary mode necessarily has to 
be in phase in 7?. = 0, as we know analytically (see Sect. 5) 
that at this point the two mode-families have strict co- and 
counter-moving character. 

Let us consider the relation between the pattern-speed 
= — — of the mode and the two rotation rates fin and 
Qp. In particular, we are interested in the region where the 
pattern-speed of the mode lies in between the rotation rates 
of the two fluids, such that it would appear prograde when 
viewed in one fluid frame and retrograde in the other. One 
can see that this "mixed" region is characterized by the con- 
dition 

K {k + mTZ) < . (95) 
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Figure 5. Avoided crossings between "ordinary" and "superfluid" (defined at 7?. = 0) inertial modes as functions of the relative rotation 
rate TZ. The common crossing-points (94) are marked by 'x' and the shaded areas corresponds to K (k + mIZ) < 0. The labels and 
K- indicate if the two fluids are in phase or in counter-phase respectively. 



This "mixed" region is indicated in Fig. 5 and Fig. 7 as a 
shaded area, and we observe that the change of the phase- 
character of modes only happens when the mode frequency 
crosses into or out of the "mixed" region. There seems to be 
no phase-change, however, if the crossing takes place via one 
of the two common crossing points (94), which are indicated 
by 'x' in these figures. We can try to understand this as 
follows: when a mode crosses into or out of the "mixed" 
region, it means that its frequency vanishes and changes sign 
in one of the two fluid frames. In general the Coriolis force 
of the corresponding fluid is nonzero in this point, therefore 
the frequency can only be zero if the fiuid ceases to move. 
The corresponding fluid eigenfunctions therefore undergo a 
sign-change, which results in the phase-change of the mode. 
In the special case where the crossing happens via one of 
the two common crossing points, however, the Coriolis-force 
does vanish at this point and subsequently changes sign, 
therefore the mode-amplitude cannot change sign and the 
crossing takes place without a phase change. 



6.3 Varying the entrainment 

In Fig. 6 we have plotted the mode-frequencies as func- 
tions of the entrainment e. Similar to the avoided crossings 



as functions of TZ shown in Fig. 5, we observe that there 
are only "pairwise" avoided crossings, i.e. between an "ordi- 
nary" and the corresponding "superfluid" mode. We further 
note that the crossing of the zero-entrainment axis (e = 0) 
is rather special, as can be understood from the discussion 
in Sect. 4.3. At e = 0, one of the two fluid-amplitudes is 
necessarily zero, and therefore the crossing of the £ = 
axis induces a phase-change between the two fluids. This is 
exactly the behaviour observed numerically for the modes 
shown in Fig. 6. 



6.4 The two-stream instability 

It was recently discovered (Andersson et al. 2002) that su- 
perfiuid systems may, quite generally, suffer a so-called "two- 
stream instability". In the present context, this instability 
would set in when the relative velocity between the two flu- 
ids exceeds a certain critical level. This mechanism was sug- 
gested as a possible mechanism for triggering pulsar glitches 
(Andersson et al. 2003). Unfortunately, the dispersion rela- 
tion for superfluid r-modes on which the analysis of Anders- 
son et al. (2003) was based is incorrect, affecting the various 
estimates for the onset and growth of the instability (for a 
detailed discussion, see Andersson et al. (2002)). 
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Figure 6. Avoided crossings between "ordinary" and "superfluid" m = 2 inertial modes as functions of entrainment e for the TV = 1 
polytropic background model and a relative rotation rate of 7^ = 0.1. The labels k;_|_ and K— refer to the phase-character of the modes. 
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respectively. The boxes indicate regions that we zoom into in Fig. 8. 
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As the general instability mechanism discussed in An- 
dersson et al. (2002) remains sound, we expect to find in- 
ertial modes that become unstable beyond a critical rela- 
tive rotation rate TZ. For the parameter values chosen for 
Fig. 5, no such instabilities were observed within the in- 
terval — 2 <TZ < 1 that was considered. However, using the 
dispersion relation (61), we can identify a more instability- 
prone region to be, for example, a proton fraction of a;p = 0.2 
and an entrainment of £ = —2. In the neutron-star core the 
entrainment e is generally expected to be positive, but a neg- 
ative entrainment is nevertheless not unphysical. Superfluid 
*He, for example, has negative entrainment, and this is also 
expected to be the case for the neutron superfiuid in the 
neutron-star crust (Carter et al. 2003). While the present 
example serves only as a consistency check and proof of 
principle, we emphasize that these parameter- values are not 
completely unphysical. In Fig. 7 we plot the frequencies of 
the lowest-order inertial modes as functions of the TZ for this 
choice of parameters. We see that now the r-mode and 
the inertial modes 2^(2) and 2l(3) do indeed undergo an in- 
stability via the merger of the "ordinary" -type mode with its 
"superfluid" counterpart. After this merger the two mode- 
frequencies are complex conjugates, which is to be expected 
from the time-symmetry of the problem. The real part of 
Ko is strictly linear in TZ in the instability region. For the 
present set of parameters, the two common crossing points 
given by (94) are (7^, kq) = (0.5,0) and (-0.6667, 1.3333). 
These points are marked by 'x' in Fig. 7. Interestingly, the 
instability point of the r-mode coincides (up to numerical 
precision ~ 10~^) with one of the common crossing points 
discussed earlier, namely the one at which the proton-fluid 
angular momentum vanishes. Using the analytic r-mode dis- 
persion relation (61), we can verify that the instability oc- 
curs exactly at the crossing point (0.5,0). However, this is 
clearly seen not to be the case for the higher order iner- 
tial modes. One might expect the instabilities to occur in 
one of the "mixed" regions, as the mode is then prograde in 



one fluid frame and retrograde in the other (cf. for example 
Pierce (1974)). This, however, is not always the case, as il- 
lustrated in Fig. 8. We see that onset of the instabilities of 
the 2-^(2) and 2^(3) modes occur slightly outside the "mixed" 
region. Given the numerical precision of < 10"^, this should 
not be due to numerical errors. This observation serves as 
a strong motivation for a study into the stability proper- 
ties of rotating multi-fiuid systems. It would be desirable to 
attempt a derivation of useful instability criteria, eg. analo- 
gous to those derived by Friedman & Schutz (1978) for the 
single fluid problem. 



7 DISCUSSION 

In this paper we have derived the equations that govern 
inertial modes of a slowly rotating superfluid neutron star 
model in the anelastic approximation. These equations are 
more general than ones used in previous studies since they 
allow for general non co-rotating backgrounds Qn 7^ fip. We 
have discussed analytically the special cases of co-rotation 
and zero entrainment. The obtained analytical results were 
then confirmed by, and thus served as important benchmark 
tests for, our numerical calculations. We studied numerically 
the dependence of the mode frequencies on the relative ro- 
tation rate and entrainment, and found avoided crossings 
between mode-pairs. The "phase character" of the modes 
was found to be rather complex, in the sense that it can 
change when crossing into or out of a "mixed region". In a 
"mixed region" the mode frequencies lie in between the two 
background rotation rates. We have also confirmed, for the 
first time in a complete mode-calculation, the existence of 
the superfluid two-stream instability. We have studied the 
onset of this instability as a function of relative rotation 
rate, and found that contrary to intuitive expectations, the 
onset can sometimes take place slightly outside the "mixed 
region" . 
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The complicated problem of oscillations of rotating 
multi-fluid systems provides many challenges that should 
inspire future work. More detailed models should allow for 
stratified backgrounds as this would be closer to a realistic 
neutron star model. Stratification is expected to lead to a 
substantially more complex character of the mode-spectrum. 
In particular, there is likely to be avoided crossings between 
all modes and the modes will no longer be of purely "in- 
phasc" or in "counter- phase" -character. This has already 
been observed in the studies by Lee & Yoshida (2003) and 
Yoshida & Lee (2003a) in the purely co-rotating case. It 
would also be interesting to move beyond both the anelastic 
approximation and the slow-rotation approximation, in or- 
der to be able to consider rapidly spinning stars. One should 
also account for the presence of an elastic crust, perhaps 
penetrated by a neutron superfluid, and include dissipative 
processes like mutual friction and bota-roactions between 
the two fluids. Another issue that needs to bo studied in 
detail is the potential gravitational-radiation instability of 
the various modes, and a suitable adaptation of the CFS in- 
stability criterion (Chandrasekhar 1970; Friedman & Schutz 
1978) to non co-rotating backgrounds. This should also help 
shed light on the two-stream instability, the true physical 
relevance of which is difficult to assess at the present time. 



ACKNOWLEDGMENTS 

We wish to thank L. Valdettaro and M. Rieutord for allowing 

us to use their LSB code, and wc further thank B. Dintrans 
and M. Rieutord for valuable discussions and help with the 
numerics. 

RP and NA acknowledge support from the EU Pro- 
gramme 'Improving the Human Research Potential and the 
Socio-Economic Knowledge Base' (Research Training Net- 
work Contract HPRN-CT-2000-00137). GLC gratefully ac- 
knowledges support from NSF grant PHYS-0140138. NA is 
grateful for generous support from a Philip Leverhulme Prize 
fellowship. 



REFERENCES 

Alpar M., Langer S., Sauls J., 1984, ApJ, 282, 533 
Andersson N., 1998, Ap.J, 502, 708 
Andersson N., Comer G., 2001, MNRAS, 328, 1129 
Andersson N., Comer G., Prix R., 2002, in preparation 
Andersson N., Comer G., Prix R., 2003, Phys. Rev. Lett. , 
90, 091101 

Batchelor G. K., 1953, Q. ,J. R. Mcteorol. Soc, 79, 224 
Bryan G., 1889, Phfl. Trans. R. Soc. Lond., 180, 187 
Carter B., 1989, in Anile A., Choquet-Bruhat M., eds. Lec- 
ture Notes in Mathematics, Vol. 1385, Rolativistic Fluid 
Dynamics (Noto, 1987). Springer- Verlag, Heidelberg, pp 
1-64 

Carter B., Chamel N., Haensel P., 2003, in preparation 
Carter B., Khalatnikov I., 1992, Ann. Phys. (N.Y.), 219, 

243 

Carter B., Langlois D., 1998, Nucl. Phys. B, 531, 478 
Chandrasekhar S., 1970, Phys. Rev. Lett. , 24, 611 
Comer G., 2002, Found. Phys., 32, 1903 
Dintrans B., Rieutord M., 2001, MNRAS, 324, 635 



Friedman J. L., Morsink S. M., 1998, Ap,J, 502, 714 
Friedman ,]. L., Schutz B. F., 1978, ApJ, 222, 281 
Langlois D., Sedrakian D. M., Carter B., 1998, MNRAS, 

297, 1189 
Lee U., 1995, A&A, 303, 515 
Lee U., Yoshida S., 2003, ApJ, 586, 403 
Lindblom L., Mendell G., 1994, ApJ, 421, 689 
Lindblom L., Mendefl G., 2000, Phys. Rev. D. , 61, 104003 
Lockitch K. H., Friedman J. L., 1999, ApJ, 521, 764 
Ogura Y., Phillips N. A., 1962, J. Atmos. Sci., 19, 173 
Pierce J., 1974, Almost All About Waves. The MIT Press 
Prix R., 2003a, submitted 
Prix R., 2003b, submitted 

Prix R., Comer G., Andersson N., 2002, A&A, 381, 178 
Prix R., Rieutord M., 2002, A&A, 393, 949 
Rieutord M., Dintrans B., 2002, MNRAS, 337, 1087 
Sedrakian A., Wasserman I., 2000, Phys. Rev. D. , 63, 

024016 

Villain L., Bonazzola S., 2002, Phys. Rev. D. , 66, 123001 
Yoshida S., Lee U., 2003a, astro-ph/0302313 
Yoshida S., Lee U., 2003b, gr-qc/0304073 



APPENDIX A: THE ANELASTIC 
APPROXIMATION 

The anelastic approximation was first introduced in atmo- 
spheric physics (Batchelor 1953; Ogura & Phillips 1962) and 
has since also been widely used in the study of stellar os- 
cillations and convection. A more detailed analysis of the 
quality and justification of this approximation in the case 
of g-modes can be found in Dintrans & Rieutord (2001) 
and Rieutord & Dintrans (2002), and it has also been used 
recently in the study of inertial modes (Villain & Bonaz- 
zola 2002). The anelastic approximation applies for modes 
with frequencies which are small compared to the inverse of 
the sound-crossing time of the star, which characterizes the 
lowest-order p-modc frequency. High frequency modes such 
as p-modes are effectively "filtered out" by this approxima- 
tion, so that only low-frequency modes like inertial modes or 
g-modes remain. We will now briefly sketch how this approx- 
imation works in the study of inertial modes of a barotropic 
star. We start from the linear perturbation equations for a 
uniformly rotating barotrope, assuming an eigenmode solu- 
tion of the form e^('^*+'^'^\ which yields 

i{uj + mQ) 5n + V -{nSv) = 0, (Al) 
i{u; + mil) Sv + 2nzxdv + V{SJl + 64>) = 0, (A2) 

We choose an average sound speed co as the natural ve- 
locity scale, and the stellar radius R as length scale, which 
implies the sound crossing time R/co as the natural time 
scale. Therefore the dimensionless mode-frequency is 

(A3) 



co/R ' 



Inertial modes have the property that their frequencies are 
of the order of n, so we introduce 

C^-= 0{1) . (A4) 

Prom the relation between pressure- and density- 
perturbations we obtain 
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Inertial modes of nan- stratified superfluid neutron stars. 



6P = p6ji = cs^ Sp , =^> nSjl = cs^ Sn . 
We write the local sound speed cs(r) as 
cs(r) = A(r)co , 



(A5) 



(A6) 



where A(r) is a function of order unity in the bulk of the 
star, but which usually vanishes at the stellar surface. As Sji 
has the dimensions of a velocity squared, the relation (A5) 
takes the following form in natural units: 



nSfj, = X (r) 5n , 



(A7) 



so Sfj. and Sn arc seen to be of the same order except close to 
the surface if A ^ 0. In this system of units, the perturbation 
equations can now be written as 

u)X~'^i{l + mC)nS'Jl + V ■ (nSv) = 0, (A8) 
u)[i{l + mOSv + 2CzxSv] + ViSjl + S4>) = 0. (A9) 

Wo restrict ourselves to modes that have low frequencies 
compared to the sound-crossing frequency cq/R, so we as- 
sume 



w < 1. 

It is straightforward to see from (A9) that 
SJ1=0{^), 

and therefore (A8) yields 
V-{nSv) = 0{Q^ A"^). 



(AlO) 



(All) 



(A12) 



In the bulk of the star, where A ~ C(l), we can therefore 
neglect the density variation Sn in the conservation equa- 
tion, leading to an error of order 0{uj^). However, in the 
boundary layer characterized by A --^ 0{io), i.e. in the region 
where the local sound-speed is of order eg co R, the er- 
ror of neglecting the compressibility of the matter becomes 
large. Nevertheless, the overall quality of the approximation 
is generally very good [see Dintrans & Rieutord (2001)], pro- 
vided this surface boundary layer is sufficiently thin, but we 
might expect the surface boundary conditions to be modi- 
fied. This is indeed the case, as (A12) now entails the regu- 
larity condition = for stellar modes with p — » at 
the surface. Therefore the surface displacement is necessar- 
ily zero in the anelastic approximation, which filters out any 
surface waves. Another consequence of this approximation 
is seen by taking the curl of (A9), which effectively elimi- 
nates the potentials 5/i and S<f> from the system of equations. 
The velocity perturbation is therefore independent of the 
pressure- and gravitational perturbation, which can both be 
determined a-posteriori from the solution and the remaining 
component of the Euler equation. The eigenmode solution is 
therefore independent of all potential perturbations, S(f>, SP 
(or equivalently S]T). Although these perturbations were not 
assumed to be zero, they are now "slaved" to the velocity 
perturbation. 



be written in the explicit form 

rWi^ + (^1+^-^) W!,- 1(1 + 1)V^=0, 

rW^ + + Wl-l{l + l)Vl=0, 

{I - l)Qi Ulr^ -{1 + 2)Q,+i Ui+^ + rn Vl 
— (1 — en)vnm'}lW^ — Eni^nmTZWl + r'tpl, 

= Ko [(1 - en)yn Wl + Enl'nlVp] , 

(/ - l)Qi U'-^ - (/ + 2)Q(+i C/^+i + mVl + rii^ 

= KO [(1 - £p)Up Wl + EpZ/p Wi] , 

{f - l)QiUt' + l{l + 2)Qi+i U'+^ 
+ {m- l{l + 1)(1 - e„)z/„m7e} Vi 

= /to + 1)(1 - e^)v^Vi -h /(/ -h l)e„i'„ Vl] , 

(/' - l)Qi U'-^ + l{l + 2)Q,+i f/^+i 
+mVl + mWl, + l{l + l)iil 
= KO [l{l + 1)(1 - ep)i^p Vl + l{l + l)£pi^p Vl] , 

{m - l{l + 1)(1 - £„)l/nr^^7^} Ul 

-h(z^ - i)Qi vt^ + i{i + 2)Q(+i yi+i 

-{l + l)QiWi-^ + lQi+iWi+^ 

= KO [1(1 + 1)(1 - £„)l/„ K + l{l + l)£nl^„ f/p] , 

mUl + {f - l)Qi Vl-^ + l{l + 2)Qi+i Vp'+i - 
-{I + l)Qi W^-^ + IQi+i Wl+^ = 
KO + ep)i^p Ul + l{l + l)£pi/p Ul] 



APPENDIX B: THE EXPLICIT OSCILLATION 
EQUATIONS 

The general system of equations (44)- (47) for the eigen- 
mode problem together with the definitions in Sect. 6 can 
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